In thk686/keitt.ssi.2019: Learning materials for Geospatial Data Analysis in R, part of the UT Summer Statistics Institute

pacman::p_load("tigris")
pacman::p_load("sf")
pacman::p_load("tidyverse")

Download the coastline data using the tigris package. This grabs data from the US Census. Here I transform to the US National Atlas coordinate system and then extract only the Atlantic coast.

coastline.sf <- coastline() %>%
  st_as_sf() %>%
  st_transform(crs = 2163)
coastline.atlantic <- coastline.sf %>%
  filter(NAME == "Atlantic")

Plot the coast. I use st_simplify to reduce the complexity of the polygons. This greatly speeds up plotting.

ggplot() +
  geom_sf(data = coastline.atlantic %>% st_simplify(dTolerance = 1e4))

Now we create a set of polygons that buffer each of the coastline line segments. We have to make the spatial set union of those polygons to make a single buffer (plus an extra for Puerto Rico).

atlantic.buffer <- coastline.atlantic %>%
  st_buffer(1e5) %>%
  st_union() %>%
  st_sf()  # units are meters

Plot them together.

ggplot() +
  geom_sf(data = atlantic.buffer %>% st_simplify(dTolerance = 1e4)) +
  geom_sf(data = coastline.atlantic %>% st_simplify(dTolerance = 1e4), color = "darkred")

Now we grab the national borders and internal divisions. It is key to transform into the same CRS.

nat <- nation() %>%
  st_as_sf() %>%
  st_transform(crs = 2163)
divs <- divisions() %>%
  st_as_sf() %>%
  st_transform(crs = 2163)

By interesecting the buffer with the boundaries, we clip off the part of the buffer that is over the ocean.

buf_int <- st_intersection(atlantic.buffer, divs)

Again plot it, this time limiting the plot area to the buffer region.

ggplot() +
  geom_sf(data = nat %>% st_simplify(dTolerance = 1e4)) +
  geom_sf(data = buf_int %>% st_simplify(dTolerance = 1e4), fill = "red", alpha = 0.2) +
  lims(x = st_bbox(buf_int)[c(1, 3)], y = st_bbox(buf_int)[c(2, 4)])

Now for all coastlines.

ggplot() +
  geom_sf(data = coastline.sf %>% st_simplify(dTolerance = 1e4))

As before, except with all coastlines.

all.buffer <- coastline.sf %>% st_buffer(1e5) %>% st_union() %>% st_sf()

ggplot() +
  geom_sf(data = all.buffer %>% st_simplify(dTolerance = 1e4)) +
  geom_sf(data = coastline.sf %>% st_simplify(dTolerance = 1e4), color = "darkred")

all_buf_int <- st_intersection(all.buffer, divs)

ggplot() +
  geom_sf(data = nat %>% st_simplify(dTolerance = 1e4)) +
  geom_sf(data = all_buf_int %>% st_simplify(dTolerance = 1e4), fill = "red", alpha = 0.2) +
  lims(x = st_bbox(all_buf_int)[c(1, 3)], y = st_bbox(all_buf_int)[c(2, 4)])

We could try to filter out all of the parts of the US not of interest using dplyr::filter.

states.sf <- states() %>% st_as_sf() %>% st_transform(crs = 2163)
# states.sf %>% filter(!NAME %in% c("Hawaii", "Puerto Rico"))

Plot the final output.

ggplot() +
  geom_sf(data = states.sf %>% st_simplify(dTolerance = 1e4)) +
  geom_sf(data = all_buf_int %>% st_simplify(dTolerance = 1e4), fill = "red", alpha = 0.2) +
  lims(x = st_bbox(states.sf)[c(1, 3)], y = st_bbox(states.sf)[c(2, 4)])

thk686/keitt.ssi.2019 documentation built on June 28, 2019, 1:28 p.m.